global google: env google
global work $google\Payments\Public
cd "$work"

***********************************
***** TABLE 7 SUMMARY STATISTICS *****
***********************************
use "$work\Stata\lab-public.dta", clear

sum rmssd rmssd_base total_ans correct_ans math confidence risk male age retailexpmths 

outreg2 using "$work\Stata-out\lab.doc", replace sum(log) wide label eqkeep(N mean sd min max) ///
		keep(rmssd rmssd_base total_ans correct_ans math confidence risk male age retailexpmths) ///
		sortvar(rmssd rmssd_base total_ans correct_ans math confidence risk male age retailexpmths) ///
		title("Laboratory experiment: Summary statistics")

/*Gender differences in HRV*/
sum rmssd_base rmssd if male==1
sum rmssd_base rmssd if male==0

********************************************************
***** FIGURE 4 HRV OVER STAGES OF THE EXPERIMENT*****
********************************************************
use "$work\Stata\lab-public.dta", clear
drop rmssd_delta
reshape long rmssd_, i(id) j(stage) string
gen stage1=stage
replace stage1="4" if stage1=="3"
replace stage1="3" if stage1=="2"
replace stage1="2" if stage1=="1"
replace stage1="1" if stage1=="math"
replace stage1="0" if stage1=="base"
destring stage1, replace
label define stage 0 "baseline" 1 "math" 2 "cashier 0-5 min" 3 "cashier 5-10 min" 4 "cashier 10-15 min"
label values stage1 stage
keep id cash rmssd_ stage1
rename rmssd_ rmssd
rename stage1 stage 

collapse (mean) meanrmssd=rmssd (sd) sdrmssd=rmssd (count) n=rmssd, by(stage cash)

*create the high and low values of the confidence interval
generate hirmssd = meanrmssd + invttail(n-1,0.025)*(sdrmssd / sqrt(n))
generate lormssd = meanrmssd - invttail(n-1,0.025)*(sdrmssd / sqrt(n))

generate stagecash = stage if cash == 1
replace stagecash = stage+7 if cash==0
twoway 	bar meanrmssd stagecash if stage==0 || ///
		bar meanrmssd stagecash if stage==1 || ///
		bar meanrmssd stagecash if stage==2 || ///
		bar meanrmssd stagecash if stage==3 || ///
		bar meanrmssd stagecash if stage==4 || ///
		rcap hirmssd lormssd stagecash, ///
		ylabel(0(10)60) ///
		ytitle("HRV") xtitle("") ///
		legend(order(1 "baseline" 2 "math task" 3 "cashier task 0-5 min " 4 "cashier task 5-10 min " 5 "cashier task 10-15 min")) ///
		xlabel( 2 "High-cash" 9 "Low-cash", noticks) ///
		plotregion(margin(l+5 r+5 b=0))
		
graph export "$work\Figures\labHRV.png" , replace as(png) width(1200) height(900)

*************************************************
**** FIGURE 5, TABLE 8 EARNINGS AND STRESS ******
*************************************************
use "$work\Stata\lab-public.dta", clear

gen earning=correct_ans*0.03
bysort cash: egen earning_mean= mean(earning) 
bysort cash: egen earning_sd  = sd(earning)  
bysort cash: gen earning_std = (earning-earning_mean)/earning_sd

label variable earning_std "earnings (standardised)"
label variable rmssd_delta "RMSSD change"
label variable cash "high-cash treatment"
label variable rmssd_base "baseline RMSSD"
label variable math "math score"
label variable confidence "confidence"
label variable risk "risk"
label variable male "male"
label variable age "age"
label variable retail_exp "retail experience"
gen excess=(total_excess>0)
label variable excess "excess change given"
gen cash_rmssd_delta=cash*rmssd_delta
label variable cash_rmssd_delta "cash x RMSSD change"

/*Figure 5*/
twoway (scatter earning_std rmssd if cash==1) (lfit earning_std rmssd if cash==1, lwidth(medthick)), ///
		ytitle("Earnings (standardized)") xtitle("HRV") ///
		title("High-cash") name("cash", replace) legend(off)
	
twoway (scatter earning_std rmssd if cash==0) (lfit earning_std rmssd if cash==0, lwidth(medthick)), ///
		ytitle("Earnings (standardized)") xtitle("HRV") ///
		title("Low-cash") name("card", replace) legend(off)
graph combine cash card, ycommon
graph export "$work\Figures\labearning.png" , replace as(png) width(1200) height(900)

/*stress difference between treatments*/
ttest rmssd_delta, by(cash)
ranksum rmssd_delta, by(cash)

ttest rmssd_delta if earning_std >0, by(cash)
ranksum rmssd_delta if earning_std >0, by(cash)

/*Table 8 regressions*/
global outopt plain coljust(lc) varlabels se starloc(1) starlevels(10 5 1) bdec(3) summstat(N) summdec(0)  summtitles(Observations)

outreg, clear(cash)
reg earning_std rmssd_delta rmssd_base math confidence risk excess male age retail_exp if cash==1, cluster(session)
outreg , $outopt merge(cash) keep(rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "High-cash")
reg earning_std rmssd_delta rmssd_base math confidence risk excess male age retail_exp if cash==0, cluster(session)
outreg , $outopt merge(cash) keep(rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "Low-cash")
reg earning_std cash rmssd_delta cash_rmssd_delta rmssd_base math confidence risk excess male age retail_exp, cluster(session)
outreg , $outopt merge(cash) keep(cash rmssd_delta cash_rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "Both treatments")
reg earning_std math confidence risk excess male age retail_exp openness conscientiousness extraversion agreeableness neuroticism, cluster(session)
outreg , $outopt merge(cash) keep(math confidence risk excess male age retail_exp openness conscientiousness extraversion agreeableness neuroticism) ctitle("VARIABLES", "Big 5")

outreg using "$work\Stata-out\lab.doc", $outopt replay(cash) title("Earning and stress") addtable replace

/*robustness - winsorize outliers*/
winsor2 earning_std, cuts(5 95) by(cash)

twoway (scatter earning_std_w rmssd if cash==1) (lfit earning_std_w rmssd if cash==1, lwidth(medthick)), ///
		ytitle("Earnings (standardized)") xtitle("HRV") ///
		title("High-cash") name("cash", replace) legend(off)
	
twoway (scatter earning_std_w rmssd if cash==0) (lfit earning_std_w rmssd if cash==0, lwidth(medthick)), ///
		ytitle("Earnings (standardized)") xtitle("HRV") ///
		title("Low-cash") name("card", replace) legend(off)
graph combine cash card, ycommon

outreg, clear(cash)
reg earning_std_w rmssd_delta rmssd_base math confidence risk excess male age retail_exp if cash==1, cluster(session)
outreg , $outopt merge(cash) keep(rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "High-cash")
reg earning_std_w rmssd_delta rmssd_base math confidence risk excess male age retail_exp if cash==0, cluster(session)
outreg , $outopt merge(cash) keep(rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "Low-cash")
reg earning_std_w cash rmssd_delta cash_rmssd_delta rmssd_base math confidence risk excess male age retail_exp, cluster(session)
outreg , $outopt merge(cash) keep(cash rmssd_delta cash_rmssd_delta rmssd_base math confidence risk excess male age retail_exp) ctitle("VARIABLES", "Both treatments")
